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DEGAS 

A  SUBROUTINE  TO  ISOLATE  MATERIAL  PACKAGES 
IN  THE  HELP  HYDRODYNAMIC  CODE 


1 .  BACKGROUND 


1,1  Use  of  the  HEW  code  at  URL 

Sophisticated  coniputer  codes  which  model  material  behaviour  are 
being  used  increasingly  as  an  adjunct  to  experimental  investigations.  "Vheir 
purpose  is  to  guide  the  experimenter  in  his  choice  of  experiments  and  to  help 
him  to  understand  the  physics  underlying  his  experiments.  Frequently  such 
computer  codes  can  be  used  to  reduce  drastically  the  number  of  difficult  or 
expensive  experiments  required  by  a  research  programme. 


At  MRL  coniputer  codes  are  being  used  to  model  the  behaviour  of 
explosives  and  explosively  driven  devices.  Computer  modelling  is  particularly 
appropriate  to  ejqilosives  research,  as  direct  measurement  of  explosive 
phenomena  requires  sophisticated  instrumentation  and  expensive  facilities. 

The  Warheads  Research  Group  at  MRL  is  using  the  US  CODE  "HELP" 
(Hydrodynamic  Elastic  Plastic)  to  assist  in  its  experimental  research  on 
metal-lined  warheads.  HELP  is  a  two-dimensional  Eulerian  hydrodynamic  code, 
originally  develoiped  for  hypervelocity  penetration  studies.  A  detailed 
description  of  the  code  may  be  found  in  Reference  [1J.  The  version  in  use  at 
MRL  (2]  was  obtained  from  the  US  Army  Ballistic  Research  Laboratory  through 
the  auspices  of  T*FCP.  BRLHELP  contains  explosive  burn  routines  and  other 
facilities  to  allcw  the  study  of  explosive  phenomena.  This,  combined  with  the 
Eulerian  geometry  which  allows  large  material  distortions  to  be  followed  by 
the  code,  made  BRLHELP  an  attractive  choice  for  the  study  of  collapsing  liner 
warheads . 
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HELP  is  being  run  on  the  CSIRONET  CDC  7600  computer,  through  the 
CS IRONODE  installed  at  MRL.  Currently  the  code  is  being  used  in  support  of  an 
experimental  investigation  of  wide-angled  shaped  charges  for  the  attack  of 
light  armour.  One  of  the  aims  of  the  experimental  programme  is  to  optimise 
parameters  such  as  the  thickness,  curvature,  and  cone  angle  of  the  shaped 
charge  liner  to  produce  coherent  projectiles  of  the  desired  shape  and 
velocity.  Hie  use  of  HELP  minimises  the  number  of  firings  required  to  obtain 
such  data.  Figure  1  shows  a  schematic  illustration  of  one  such  device,  and 
Figure  2  shows  an  experimental  flash  radiograph  of  the  ensuing  projectile  at 
200  microseconds. 


1.2  Common  Difficulties  Encountered  in  Running  the  Help  Code 

Obtaining  the  best  performance  from  a  sophisticated  corrputer  code 
such  as  HELP  is  as  much  art  as  science.  Only  experience  can  tell  the  user 
what  values  of  input  parameters  will  result  in  the  most  stable  calculations 
for  a  given  problem  type.  However  some  difficulties  with  HELP  are  more 
intrinsic,  and  can  be  anticipated  with  most  problems. 

One  of  the  most  frequent  causes  of  instability  in  the  HELP  code  is 
the  pressure  iteration.  Hie  pressure  in  a  single-material  cell  is  computed 
directly  from  the  material  density  and  internal  energy  using  the  equation  of 
state.  Hie  Eulerian  logic  of  the  code,  however,  allows  multi-material  (mixed) 
cells,  and  the  pressures  in  these  cells  are  derived  by  iteration.  Hie 
densities  of  the  component  materials  of  a  mixed  cell  are  varied,  subject  to 
filling  the  cell  exactly,  until  the  conponent  pressures  converge  to  a  common 
value,  taken  as  the  cell  pressure.  Hiis  iteration  can  lead  to  unrealistic 
answers  when  materials  in  a  mixed  cell  have  very  disparate  densities  and/or 
masses. 

When  the  pressure  iteration  fails  under  these  conditions  (even 
though  convergence  is  obtained)  a  "rogue”  cell  is  created  where,  usually,  a 
small  amount  of  material  has  physically  unlikely  properties.  Unfortunately 
although  one  such  rogue  cell  may  be  of  little  concern  to  the  user,  it  is 
likely  to  affect  the  overall  functioning  of  the  code,  for  example  via  the 
time-step  calculation.  Hie  normal  procedure  in  such  a  case  is  surgery.  Hie 
code  senses  trouble  and  steps,  relevant  quantities  in  the  rogue  cell  are 
adjusted  by  the  user,  and  the  code  is  restarted  from  the  binary  dunp  file. 

Hiis  type  of  procedure  is  common  for  adjusting  minor  hiccups  in  the  running  of 
large  hydrodynamic  codes. 

Experience  with  the  use  of  HELP  to  model  explosive-metal 
interactions  has  shown  that  the  above  problems  occur  quite  frequently  when  a 
mixed  cell  contains  low  density  gaseous  detonation  proAicts  as  well  as 
condensed  material.  Late  in  a  calculation,  both  the  pressures  and  densities 
of  these  gases  are  decreasing  with  time,  exacerbating  the  problem.  Here  the 
approach  of  fixing  a  rogue  cell  and  continuing  is  inappropriate,  as  the 
condition  will  recur.  A  better  approach,  assuming  the  gases  are  no  longer 
making  physically  significant  contributions  to  the  problem,  is  to  remove  them 
altogether.  (For  most  explosives  problems  1  kPa  would  be  a  reasonable 


pressure  threshold. )  This  approach  has  the  additional  benefit  of  increasing 
the  run-time  economy  of  the  code. 


2.  DEGAS  -  DESCRIPTION  AND  USE 


The  purpose  of  the  subroutine  "DEGAS"  is  to  remove  explosives  gases 
from  the  grid  of  a  HELP  calculation.  The  logic  of  HELP  is  such  that  it  is  not 
particularly  easy  to  "fool"  the  code  into  believing  that  material  has 
disappeared  (for  example  by  zeroing  its  equation  of  state).  The  main  reasons 
for  this  are  the  interface  cells,  which  have  their  own  quite  complicated 
logic,  and  the  massless  tracer  particles,  which  are  sensed  by  the  code  as 
material  boundaries. 


For  simplicity,  DEGAS  does  not  remove  only  gases,  but  rather 
isolates  material  package  one,  removing  all  other  materials  and  associated 
array  quantities  from  the  calculational  grid.  This  will  almost  always  satisfy 
those  using  the  code  for  warhead  studies,  where  the  projectile  is  of  prime 
interest,  and  other  materials  such  as  the  expanded  explosive,  as  well  as  any 
casings,  wave-shaping  inerts  etc,  may  be  discarded  once  their  influence  on  the 
projectile  has  ceased.  DEGAS  can  however  be  modified  by  the  knowledgeable 
user  to  remove  selective  material  packages  as  required. 

Appendix  1  contains  an  UPDATE  run  [3]  to  add  DEGAS  to  the  HELP  code, 
including  a  Fortran  listing.  It  also  includes  a  minor  correction  to  the  code, 
IDENT  MOVERR,  which  prevents  the  tracers  at  the  end-points  of  material 
packages  from  being  stranded.  (The  use  of  DEGAS  made  this  error  more  likely 
to  manifest  itself.) 


3.  EXAMPLES  OF  THE  USE  OF  DEGAS 


Reference  (2)  details  a  sample  problem  for  use  with  HELP,  involving 
a  hemispherical  shaped  charge,  1.5  inch  diameter.  This  sample  problem  was 
used  during  the  installation  of  HELP  at  MRL,  and  Figure  3  shows  outlines  of 
the  copper  liner  towards  the  end  of  one  such  run.  Wiis  run  terminated  after 
approximately  17.5  microseconds,  due  to  problems  associated  with  the  pressure 
iteration  in  a  mixed  cell.  BRL  appear  to  have  had  similar  experiences  as  the 
output  in  [2]  stops  at  18  microseconds.  When  the  same  problem  was  restarted 
using  DEGAS,  as  would  be  expected  the  code  ran  for  much  longer.  Figure  4 
gives  output  for  up  to  24  microseconds,  and  shows  that  stretching  of  the  liner 
has  increased  significantly.  (In  fact  the  problem  ran  to  over  60 
microseconds,  but  results  after  24  microseconds  were  of  little  value  as  (i) 
the  liner  had  virtually  stopped  stretching  and  (ii)  the  tracers  started  to 
generate  physically  silly  results  near  the  base  of  the  liner. ) 
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Figure  5  shews  tracer  outlines  from  a  HELP  run  on  the  problem  shewn 
in  Figure  1.  Again  this  run  stopped  because  of  problems  in  a  mixed  cell 
containing  metal  and  expanded  HE,  and  again  this  occurred  before  the  metal 
liner  had  reached  a  stable  shape.  With  the  geometry  used  (4  inch  diameter 
charge)  both  experiment  and  theory  predicted  that  the  liner  would  stabilise 
somewhere  between  150  and  200  microseconds.  Figure  6  shows  some  of  the  output 
from  a  restart  of  this  problem  using  DEGAS.  The  liner  reaches  a  stable  shape 
at  about  160  microseconds. 

The  "tails"  on  the  metal  liner  derive  from  failed  material,  and  are 
no  longer  behaving  in  a  physically  realistic  manner.  The  manual  removal  of 
this  material  from  the  code  is  quite  complex,  so  a  subroutine  CLIP  has  been 
written  to  do  this.  CLIP  is  quite  similar  to  DEGAS  in  that  most  of  the 
subroutine  is  taken  up  with  "housekeeping"  so  that  the  code's  accounting 
systems  are  satisfied.  Appendix  2  contains  an  UPDATE  run  [3]  to  add  CLIP  to 
the  HELP  code,  and  includes  a  Fortran  listing. 

Figure  7  shows  some  output  from  a  restart  of  the  above  problem, 
using  DEGAS  and  CLIP.  Cropping  the  metal  liner  in  this  case  produces  output 
that  is  more  readily  understood:  it  also  decreases  the  code’s  runtime  and 
prevents  tracers  from  crossing  each  other.  To  some  extent  it  will  also  lead 
to  more  reliable  calculations,  as  rogue  cells  from  failed  (and  physically 
unrealistic)  portions  of  the  material  are  prevented  from  controlling  the 
timestep. 
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COPPER  LINER 


All  dimensions  in  inches 


FIGURE  1.  Schematic  of  Wide-angled  Shaped  Charge. 
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FIGURE  3«  Hemispherical  copper  liner  run. 
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FIGURE  6.  (contd.  ) 
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FIGURE  7.  Wide-angled  Shaped  Charge  Restart  with  DEGAS  and  CLIP. 


APPENDIX  1 


AN  UPDATE  RUN  TO  ADD  DEGAS  TO  THE  HELP  CODE 


MAKE6, PI  000,  T1 00. 

GETSET(SN=DFC4405) 

ATTACH(OLDPL,MRLHLP5,SN-DFC4405,ID=DFCDLS) 

REQUEST! NEWPL,*PF,SN=DFC4405) 

UPDATE! F,L=A1 234, N) 

REQUEST! LGO ,  •  PF , SN=DFC440  5 ) 

FTN! I=C0MPILE,LCM=1 ) 

CATALOG!  NEWPL,$MRLHLP6  FOR  EXPORT? , SN=DFC4405, ID=DFCDLS,RP=1 00) 

CATALOG! LGO, ?NEW  LGO  FOR  MRLHLP6? , SN=DFC4405 , ID=DFCDLS ,RP=1 00) 

•IDENT  MOVERR 
♦DELETE  MOVTCR .163,164 

C  THIS  FIX  BY  DLS  TO  PREVENT  TRACERS  1  AND  NN  GETTING  STUCK. 

C  IT  APPEARS  THIS  WAS  INTENDED  ORIGINALLY,  AT  LEAST  FOR  THE 

C  L=NN  CASE.  THE  PROBLEM  BECAME  MORE  LIKELY  TO  ARISE  WITH 
C  THE  USE  OF  THE  DEGAS  ROUTINE  TO  ISOLATE  PACKAGE  ONE. 

C 

IF(L.EQ.NN)NT=-1 

IF! TX!N , L+NT) .LE. -1 000. )GO  TO  1505 
C  *  *  *NOTE  PREVIOUSLY  THE  PROG  SKIPPED  TO  1505  IF  L  WAS  1  OR  NN. 

C  ***ALSO  TX  FOR  L+1  AND  L-1  WERE  TESTED,  WITHOUT  USING  NT=-1 . 

C 

•IDENT  NOGAS 
♦DELETE  FIT. 23 

1  ICSTOP  ,PIDY  , TOPMU  , RTMU  , NOGAS  ,NUMEREZ  , ETH 

•DELETE  INPUT. 35 

1  SIEMIN,  STAB,  TSTOP,  NOGAS 

•INSERT  INPUT. 90 

C 

C  CALL  DEGAS  ON  RESTART  IF  NOGAS  EQUALS  1 

C 

IF  ( PK( 3 ).LT.O.. AND. NOGAS. EQ.1) CALL  DEGAS 
•ADDFILE 


•DECK  DEGAS 

SUBROUTINE  DEGAS 


C 

C 

C 

C 

C 

C 

C 

C 

C 

c 

c 

C 

c 

c 


*  THIS  SUBROUTINE  ISOLATES  PACKAGE  1  ON  A  RESTART  RUN,  AND 

*  RESTORES  THE  PROGRAM  LOGIC  TO  ALLOW  CALCULATION  TO  PROCEED. 

*  DEGAS  IS  CALLED  WHEN  THE  Z- BLOCK  VARIABLE  NOGAS  (ZC 1 1 ) )  IS  1, 

*  FOLLOWING  A  RESTART.  DEGAS  IS  USED  TO  REMOVE  DETONATION 

*  PRODUCTS  WHICH  MAY  BE  ADVERSELY  AFFECTING  PRESSURE  ITERATIONS 

*  LONG  AFTER  THEIR  CONTRIBUTION  TO  THE  PROBLEM  HAS  CEASED. 

* 

*  DEGAS  ASSUMES  THAT  THE  PACKAGE  YOU  WISH  TO  ISOLATE  (USUALLY  A 

*  WARHEAD  LINER  OR  SIMILAR)  IS  PACKAGE  ONE. 

* 

*  AUTHOR  DAVID  L.  SMITH,  AUGUST  1982. 


A1-2 


C  ************************* 

c 

♦CALL  COMDK 
C 

WRITE(6,600) 


C 

Q  ************************* 

C  *  LOOP  TO  ADJUST  THEORETICAL  ENERGY  FOR  REMOVAL  OF  MATERIALS 


Q  ************************* 

c 

ESUM=0. 

DO  30  K=2 ,  KMAX 
M=IABS(MFLAG(K) ) 

IF(M.EQ.1 )GOTO  30 
IF(M.LT. 1 00) GOTO  20 
M=M- 1 00 

DO  10  N=2,NMAT 

ESUM=ESUM+XMASS (N, M) * ( SIE(N, M) +TKEGM ( N,  M) ) 

1 0  CONTINUE 
GO  TO  30 

20  ESUM=ESUM+AMX(K)*(TKEG(L)+AIX(K) ) 

30  CONTINUE 

0  ************************* 

C  *  NOW  DEDUCT  ESUM  FROM  ETH,  THE  THEORETICAL  ENERGY  OF  THE  GRID.  * 

C  ************************* 

C 

ETH= ETH- ESUM 

0  ************************* 

C  *  BEGIN  K-LOOP  FOR  PURE  CELLS  ONLY  * 

Q  ************************* 

C 

DO  100  K= 1 , KMAX 
MFK=MFLAG(K) 

IF(MFK.GT. 1 00) GOTO  100 
IF(MFK.EQ. 1 )GOTO  100 

C  ************************* 


C  *  ZERO  ALL  PURE  CELLS  EXCEPT  PACKAGE  1  * 

q  ************************* 


C 

U(K)=0. 

V(K)=0. 

AMX(K)=0. 

AIX(K)=0. 

TKEG ( K ) =0 . 

P(K)=0. 

MFLAG(K)*0 

DETIM(K)=0. 

STRSZZ(K) =0. 
STRSRR(K)=0. 

STRSRZ(K) =0. 

C  END  LOOP  FOR  PURE  CELLS 
1 00  CONTINUE 
C 


Al-3 


Q  ************************* 

c  *  BEGIN  LOOP  FOR  MIXED  CELLS  * 

C  ************************* 

C 

DO  200  K=1 ,KMAX 
MFK=MFLAG(K) 

I F ( MFK . LE . 1 00 ) GOTO  200 
M=MFK-1 00 

Q  ************************* 

c  *  ZERO  ALL  MATERIALS  IN  MIXED  CELLS  EXCEPT  PACKAGE  1  * 

C  ************************* 

c 

DO  210  N=2 , NMAT 
XMASS(N,M)=0. 

RHO(N,M)=0. 

SIE(N, M) =0. 

TKEGM  ( N ,  M  )  =0 . 

US(N,M)=0. 

VS(NfM)=0. 

SAMPY(N,M)=0. 

SGAMC(N,M)=0. 

SAMMY(N,M)=0. 

SAMMP(N,M)=0. 

210  CONTINUE 

Q  ************************* 

C  *  CHECK  FOR  PRESENCE  OF  PACKAGE  1  * 

Q  ************************* 

C 

IF(RHO(1 ,M).GT.O. )GOTO  250 
U(K)=0. 

V(K)=0. 

AMX(K)=0. 

AIX(K)=0. 

TKEG(K)=0. 

P(K)=0. 

STRSZZ(K)=0. 

STRSRRC  K) =0. 

STRSRZ(K)=0. 

DETIM  ( K)  =0 . 

c  ************************* 

C  *  SET  FLAG  TO  RELEASE  MIXED  CELL  * 

C  *******************  ‘  *  *  *  *  * 

C 

RHO(1 ,M)=-1 . 

c  ************************* 

C  *  CONVERT  TO  PURE  VOID  CELL  * 

C  ***.************»**•*•*** 

C 

MFLAG(K)»0 

C 

GOTO  200 
C 

250  CONTINUE 


i 


A1  -4 


C  ************************* 

C  *  CELL  IS  KNOWN  TO  CONTAIN  PACKAGE  1  * 

C  *  * 

C  *  CELL  MUST  ALSO  NOW  CONTAIN  VOID  PACKAGE  * 

C  *  HENCE  FLAG  VOID  AS  1,  USING  NEW  NVOID(2).  * 

C  ************************* 

c 

RHO(2,M)=1 . 

Q  ************************* 

C  *  CHANGE  CELL  PROPS  TO  PACKAGE  1  PROPS  * 

Q  ************************* 

c 

AMX(K)=XMASS(1 ,M) 

AIX(K)=SIE(1 ,M) 

U(K)=US( 1 ,M) 

V(K)=VS(1  ,M) 

TKEG ( K ) =TKEGM ( 1 ,M) 

C  END  LOOP  TOR  MIXED  CELLS 

200  CONTINUE 

C  ************************* 

C  *  REMOVE  SLIPLINE(S)  FROM  PROBLEM  * 

Q  ************************* 

c 

DO  300  L=1  ,NMXCLS 
300  THETA(L)=-1 .0 
C 

DO  310  N=1,NMAT 
MASTRD ( N ) =0 
310  NSLAVD(N)=0 

Q  ************************* 

C  *  ADJUST  TRACERS.  * 

C  *  * 

C  *  ( 1 )  REMOVE  TRACERS  FROM  OLD  PACKAGES  * 

C  ************************* 

C 

DO  400  N=2,NVOID 
DO  400  L=1  ,NTPMX 
TX(N,L)=0. 

TY(N,L)=0. 

400  CONTINUE 

C  ************************* 

C  *  (2)  PUT  TRACERS  BACK  FOR  VOID  (REFLECT  PACKAGE  1)  * 

Q  ************************* 

c 

NNN=NMP ( 1  ) - 1 
DO  500  N=1 ,NNN 
NN2»NNN-N+1 
TX(2,N)»TX(1 ,NN2) 

TY(2,N)-TY(1 ,NN2) 

500  CONTINUE 

TX(2,NNN+1 )— 1000. 

TY(2,NNN+1)«0. 

NMP ( 2 ) “NMP ( 1  ) 


on  oooo  on 


Al-5 


c  ************************* 

C  *  ADJUST  Z-BLOCK  VARIABLES  * 

C  ************************* 

c 

NMAT=1 

NOSLIP=1 

NSLD=0 

NTCC=0 

ICLADD=-1 

Q  ************************* 

C  *  ADJUST  BLANK  COMMON  VARIABLES  * 

************************* 

MOS=0 
NVOID=2 
RM0M=0. 

ZMOM=0. 

************************* 

*  RESET  NOGAS  IN  Z-BLOCK  WHEN  DEGAS  COMPLETE  * 

************************* 

N0GAS=0. 


WRITE( 6,610) 

C 

RETURN 

C 

600  FORMAT(1H1 ,40H**»  SUBROUTINE  DEGAS  HAS  BEEN  CALLED  ***) 
610  FORMAT( 1  HO, 33H***  SUBROUTINE  DEGAS  COMPLETE  *•*) 

C 


END 


APPENDIX  2 


AN  UPDATE  RUN  TO  ADD  CLIP  TO  THE  HELP  CODE 

COMPIL ,  PI  00 ,  T1 00. 

GETSET{ SN=DFC4405) 

ATTACH!  OLDPL, $NEW  MRLHLP8S ,SN=DFC4405,ID=DFCDLS) 

REQUEST! NEWPL, *PF , SN=DFC4405) 

REQUEST! LGO, *PF,SN=DFC4405) 

UPDATE (F,N,L=A1  234) 

FTN( LCM=1 , 1 -COMPILE) 

CATALOG! NEWPL, $MRLHLP9$ .SN-DFC4405, ID-DFCDLS , RP= 900) 

CATALOG! LGO, $LGO  FOR  MRLHLP95 , ID-DFCDLS , SN-DFC4405 ,RP=900) 

♦IDENT  CLIPIT 
♦INSERT  INPUT. 91 
C 

C  CALL  CLIP  ON  RESTART  IF  NIMAX  GT  ZERO 
C 

IF(PK(3) .LT.O. . AND. NIMAX , GT. 0 ) CALL  CLIP 

C 

C 

♦DELETE  FIT.  25 

3  IOl  , SLPNDX  , SLPNDY  ,DMIN  .NIMAX  ,DTNA  ,CVIS 

•DELETE  INPUT. 33 

9  NVRTEX ,  NIMAX,  PLGOPT,  PIMMIN,  PMIN,  PRCNT, 

♦ADDFILE 
♦DECK  CLIP 

SUBROUTINE  CLIP 

C  ************************* 

C  *  ROUTINE  TO  CLIP  EDGES  OF  PACKAGE  1 ,  BY  LIMITING  IN  X-DIRECTION  * 

C  *************•***•»***••• 

c 

♦CALL  COMDK 

XSTOP- FLOAT! NIMAX ) 

KOUNT-O 

C  ****»**********»*•*•**••• 
C  *  LOOP  TO  FIND  STRADDLING  PAIRS  OF  TRACERS  * 

C  ‘  *  ‘  *******  ‘  *  *  ‘  •****»»•**  * 

c 

DO  10  I=»2,NTPMX 
TX 1 =TX ( 1 ,1-1 ) 

TX2=TX(1 ,1) 

TY1=TY(1 ,1-1) 

TY2»TY(1  ,1) 

IF(TX1 . LE. XSTOP . AND. TX  2 . LE . XSTOP ) GOTO  10 
IF(TX1 , GE. XSTOP. AND . TX2 .GE . XSTOP )GOTO  10 
C  ************************* 
C  *  HHNCE  THIS  IS  A  STRADDLE  CASE  * 

C  ************************* 

c 

NBIGT-1 

IF(TX1 .GT.TX2)NBIGT»I-1 


C  * 
C  * 
C  * 


INTERPOLATE  Y  VALUE  FOR  NEW  TRACER 


o  n 


pa-2 


c 

IF( IVARDX. EQ. 0. AND. IVARDY.  EQ. 0 ) GOTO  1 1 

C  ************************* 

C  *  CONVERT  TO  CM  UNITS  IF  NECESSARY.  * 

C  ************************* 

C 

TX 1  — XCTOP ( TX 1 ) 

TX2=XCTOP(TX2) 

TY1=YCTOP(  TY1 ) 

TY2=YCTOP(TY2) 

XSTOPC=XCTOP ( XSTOP ) 

TYNEW=TY1 +( XSTOPC-TX1 ) * ( TY2-TY1 )/{TX2-TX1 ) 

TY( 1 , NBIGT ) =YPTOC ( TYNEW ) 

TX( 1 , NBIGT )=XSTOP 
GOTO  10 

11  TY( 1 , NBIGT ) =TY 1 + ( XSTOP-TX 1 )*(TY2-TY1 )/(TX2-TX1 ) 

TX( 1 , NBIGT )=XST0P 
1 0  CONTINUE 

NNN=NMP(  1  ) 

DO  100  1=1 ,NNN 

IF( TX( 1 , I ) .GY. XSTOP) GOTO  50 

KOUNT=KOUNTf  1 

TX (  3  , KOUNT ) -TX ( 1 ,1) 

T  Y  (  3 , KOUNT ) =T Y (1,1) 


50 

TX ( 1 ,I)=0. 

TY(1 , I ) =0. 

100 

CONTINUE 

C 

***** 

****** 

*  *  *  * 

* 

*  *  * 

C 

* 

NOW  TRANSFER 

TRACERS  BACK 

TO 

PACKAGE 

c 

***** 

****** 

*  *  *  * 

* 

*  *  * 

C 

DO  110  1=1, KOUNT 
TX(1 , I)=TX( 3,1) 

TY(  1  , 1 )  =T  Y  (  3 , 1 ) 

1 1 0  CONTINUE 
C 

NMP( 1 )=KOUNT 

Q  ************************* 

C  *  TRACERS  FOR  VOID  . . .  REFLECT  PACKAGE  1  * 

************************* 


NNN=NMP( 1 ) -1 
DO  490  1=1  , NTPMX 
TX  (  2 , 1 )  =0 . 

490  TY( 2 , I )=0. 

DO  500  N= 1 ,NNN 
NN2=NNN-N+1 
TX ( 2 , N ) =TX ( 1 ,NN2 ) 
TY( 2 ,N)=TY( 1 ,NN2) 
500  CONTINUE 

TX( 2 ,NNN+1 ) =-1000. 
TY ( 2 , NNN+ 1 ) =0 . 
NMP(2)=NMP( 1 ) 
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C  ************************* 

C  *  ZERO  ALL  CELLS  PAST  NIMAX  * 

C  ************************* 


C 

ESUM=0. 

III=NIMAX+1 
DO  200  J=1 ,  JMAX 
DO  200  I=III , IMAX 
K=( J-1 )*IMAX+I+1 

C  ************************* 

C  •  ADJUST  THEORETICAL  ENERGY  FOR  REMOVAL  OF  MATERIALS  * 

C  ************************* 

C 

M-IABS(MFLAG(K) ) 

IF(M.LT.100)GOTO  230 
M=M- 100 

DO  210  N=1 ,NMAT 

ESUM=ESUM+XMASS(N,M)*(SIE(N,M)+TKEGM(N,M) ) 

210  CONTINUE 
GOTO  200 

230  ESUM=ESUM+AMX(K)*(TKEG(K)+AIX(K) ) 

200  CONTINUE 


c 

* 

* 

*  *  * 

* 

* 

* 

* 

*  *  * 

*  *  * 

* 

*  * 

*  * 

* 

* 

* 

* 

* 

c 

* 

DEDUCT 

ENERGY  FROM 

THEORETICAL 

TOTAL, 

ETH 

* 

c 

* 

* 

*  *  * 

* 

* 

* 

* 

*  *  * 

*  *  W 

* 

*  * 

*  * 

* 

* 

* 

• 

* 

c 

ETH= 

:ETH-ESUM 

c 

* 

* 

*  *  * 

* 

* 

* 

* 

*  *  * 

*  *  * 

* 

*  * 

*  * 

* 

* 

* 

• 

* 

c 

* 

ZERO  ALL  CELLS 

# 

c 

* 

* 

*  *  * 

* 

* 

* 

* 

*  *  * 

*  *  * 

* 

*  * 

*  * 

* 

* 

* 

* 

* 

C 

DO  250  J=1,JMAX 
DO  250  I— III , IMAX 
K=( J-1 )*IMAX+I+1 
U(K)=0. 

V(K)=0. 

AMX(K)=0. 

AIX(K)=0. 

TKEG(K)=0. 

P(K)=0. 

STRSZZ(K)=0. 

STRSRR(  K)  “0 . 
STRSRZ(K)=0. 
MFK=MFLAG(K) 

MFLAG(K)«0 

IF(MFK.GT. 100) GOTO  260 
GO  TO  250 
260  M-MFK-100 

DO  270  N-1 ,NMAT 

XMASS(N,M)»0. 

RHO{N,M)-0. 

SIE(N,M)-0. 

TKEGM(N,M)«0. 


'f 


% 


A2-4 


US(N,M)=0. 

VS(N,M)=0. 

SAMP¥(N,M)=0. 

SGAMC(N,M)=0. 

SAMMY(N,M)=0. 

SAMMP(N,M)=0. 

270  RHO(1,M)=-1. 

250  CONTINUE 

C  ************************* 

C  *  REDUCE  ACTIVE  GRID  * 

Q  ************************* 

c 

I1-NIMAX+2 

c  ************************* 
C  *  RESET  NIMAX  * 

C  ************************* 

c 

NIMAX=0 
CALL  ADDTCR 

C 

RETURN 

END 
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